library(foreign)
library(ggeffects)

#laptop
setwd("/Users/joelsievert/Dropbox/Research/Legislative Responsiveness/National Road/data/PRQ Replication")



h21 <- read.csv("h21_extend.csv")
h23 <- read.csv("h23_extend.csv")


#######################################
#figure 1a: distance to national road#
######################################

#washington road
m1 <- glm(rc_natl ~ jackson  + dist_road + dist_natl + slave_p +  distrib_cmte, family = binomial, data = h21)

##23rd
m5 <- glm(rc ~ jackson   + dist_road + dist_natl + slave_p + distrib_cmte , family = binomial, data = h23)



p1 <- ggpredict(m1, terms = "dist_road[1, 1.2, 1.4, 1.6, 1.8, 2, 2.2, 2.4, 2.6, 2.8, 3]", condition = c(jackson = 1, dist_natl = 3.3, slave_p = 0, distrib_cmte = 0), ci_level = 0.9)
p2 <- ggpredict(m5, terms = "dist_road[1, 1.2, 1.4, 1.6, 1.8, 2, 2.2, 2.4, 2.6, 2.8, 3]", condition = c(Jackson = 1, dist_natl = 9.6, slave_p = 0, distrib_cmte = 0), ci_level = 0.9)

par(mfrow = c(2,2))

#Washington Road
par(mar = c(5, 4, 3, 2), bg = "white")
plot(p1$x, p1$predicted, ylim = c(0, 1), xlim = c(1, 3), type = 'n', xlab='', ylab='', axes= FALSE )

abline(v = seq(1, 3, .1), col = 'grey85', lwd = 80)
abline(v = seq(1, 3, .2), h = seq(0, 1, 0.1), col = 'white')
lines(p1$x, p1$predicted, lwd = 3)
lines(p1$x, p1$conf.low, lwd = 3, lty = 2)
lines(p1$x, p1$conf.high, lwd = 3, lty = 2)

#code puts axis text and an angle
axis(side = 2, at = seq(0, 1, 0.1), labels = seq(0, 1, 0.1), las = 2, lwd = 0 ) 
axis(side = 1, at = seq(1, 3, 0.2), lwd = 0)

mtext("Distance to Nat'l Road (100s miles)", side = 1, line = 2.75)
mtext("Predicted Probability of Support", side = 2, line = 2.5)

mtext("Washington Road", side = 3, line = 0.5, cex = 1, font = 2)

#New Orleans Route
par(mar = c(5, 4, 3, 2), bg = "white")
plot(p2$x, p2$predicted, ylim = c(0, 1), xlim = c(1, 3), type = 'n', xlab='', ylab='', axes= FALSE )

abline(v = seq(1, 3, .1), col = 'grey85', lwd = 80)
abline(v = seq(1, 3, .2), h = seq(0, 1, 0.1), col = 'white')
lines(p2$x, p2$predicted, lwd = 3)
lines(p2$x, p2$conf.low, lwd = 3, lty = 2)
lines(p2$x, p2$conf.high, lwd = 3, lty = 2)

#code puts axis text and an angle
axis(side = 2, at = seq(0, 1, 0.1), labels = seq(0, 1, 0.1), las = 2, lwd = 0 ) 
axis(side = 1, at = seq(1, 3, 0.2), lwd = 0)

mtext("Distance to Nat'l Road (100s miles)", side = 1, line = 2.75)
#mtext("Predicted Probability", side = 2, line = 2.75)

mtext("New Orleans Route", side = 3, line = 0.5, cex = 1, font = 2)



#######################################
#figure 1b: distance to extension#
######################################

#washington road
m3 <- glm(rc_mays ~ jackson + dist_road + dist_mays + slave_p +  distrib_cmte, family = binomial, data = h21)

##23rd
m5 <- glm(rc ~ jackson   + dist_road + dist_natl + slave_p + distrib_cmte , family = binomial, data = h23)



p3 <- ggpredict(m3, terms = "dist_mays[3, 3.5, 4, 4.5, 5, 5.5, 6, 6.5, 7, 7.2, 7.5]", condition = c(jackson = 1,  dist_road = 2.1, slave_p = 0, distrib_cmte = 0), ci_level = 0.9)
p4 <- ggpredict(m5, terms = "dist_natl[6, 6.9, 7, 8, 9, 10, 11, 12]", condition = c(jackson = 1, dist_road = 2.4, slave_p = 0, distrib_cmte = 0), ci_level = 0.9)


#Maysville Road
par(mar = c(5, 4, 3, 2), bg = "white")
plot(p3$x, p3$predicted, ylim = c(0, 1), xlim = c(3, 7.5), type = 'n', xlab='', ylab='', axes= FALSE )

abline(v = seq(3, 7.5, .1), col = 'grey85', lwd = 80)
abline(v = seq(3, 7.5, 0.5), h = seq(0, 1, 0.1), col = 'white')
lines(p3$x, p3$predicted, lwd = 3)
lines(p3$x, p3$conf.low, lwd = 3, lty = 2)
lines(p3$x, p3$conf.high, lwd = 3, lty = 2)

#code puts axis text and an angle
axis(side = 2, at = seq(0, 1, 0.1), labels = seq(0, 1, 0.1), las = 2, lwd = 0 ) 
axis(side = 1, at = seq(3, 7.5, 0.5), lwd = 0)

mtext("Distance to Extension (100s miles)", side = 1, line = 2.75)
mtext("Predicted Probability of Support", side = 2, line = 2.5)

mtext("Maysville Road", side = 3, line = 0.5, cex = 1, font = 2)

#New Orleans Route
par(mar = c(5, 4, 3, 2), bg = "white")
plot(p4$x, p4$predicted, ylim = c(0, 1), xlim = c(6, 12), type = 'n', xlab='', ylab='', axes= FALSE )

abline(v = seq(6, 12, .1), col = 'grey85', lwd = 80)
abline(v = seq(6, 12, 1), h = seq(0, 1, 0.1), col = 'white')
lines(p4$x, p4$predicted, lwd = 3)
lines(p4$x, p4$conf.low, lwd = 3, lty = 2)
lines(p4$x, p4$conf.high, lwd = 3, lty = 2)

#code puts axis text and an angle
axis(side = 2, at = seq(0, 1, 0.1), labels = seq(0, 1, 0.1), las = 2, lwd = 0 ) 
axis(side = 1, at = seq(6, 12, 1), lwd = 0)

mtext("Distance to Extension (100s miles)", side = 1, line = 2.75)
#mtext("Predicted Probability", side = 2, line = 2.75)

mtext("New Orleans Route", side = 3, line = 0.5, cex = 1, font = 2)

